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We present a quantum-number projection technique wliich enables us to exactly treat spin, mo- 
mentum and other symmetries embedded in the Hubbard model. By combining this projection 
technique, we extend the path-integral renormahzation group method to improve the efficiency of 
numerical computations. By taking numerical calculations for the standard Hubbard model and 
the Hubbard model with next nearest neighbor transfer, we show that the present extended method 
can extremely enhance numerical accuracy and that it can handle excited states, in addition to the 
ground state. 



PACS numbers: 71.10.Fd,71.10.-w,02.70.-c,71.15.Qe 

I. INTRODUCTION 

Quantum many-body systems often possess several 
symmetries. For example, the Hubbard model preserves 
total spin, total momentum, and some geometrical sym- 
metries on a lattice. It is crucially important to iden- 
tify the symmetry and quantum numbers in understand- 
ing the nature of the ground state, where a symmetry 
breaking, for example, often occurs in the thermody- 
namic limit. The symmetry should be restored in finite 
size systems. However even in finite-size systems, the 
ground state and excitation spectra reflect the natures in 
their thermodynamic limits. Their excitation spectra and 
spectroscopic properties are resulted from eigenstates of 
specified quantum numbers and play crucial roles in eluci- 
dating the nature of low-energy phenomena in condensed 
matter physics. 

To investigate quantum many-body problems, quan- 
tum Monte Carlo (QMC) approaches have been one of 
useful methods and can give ground state properties if 
there is no minus sign problem (jy. However, they can 
not fully take an advantage of symmetry explicitly and 
excitation spectra have not been well explored. Further- 
more, if the minus sign problem becomes serious as in 
the case of the Hubbard model on a non-bipartite lat- 
tice, QMC methods do not give a well convergent result. 
Although the exact diagonalization methods handle the 
whole excitation spectra, tractable system size is severely 
limited. 

In nuclear structure physics, symmetry plays a primar- 
ily important role. For instance, as nucleus is a finite sys- 
tem, rotational symmetry is specially important. There- 
fore symmetry has been continuously focused in order 
to solve nuclear quantum many-body problems. There 
are several ways to handle symmetries in nuclear struc- 
ture problems. Among them, the projection technique 
is powerful in the respect to broken symmetry and its 
restoration. In nucleus, a mean field solution is con- 
sidered as the first approximation but it violates most 
of symmetries, i.e., total angular momentum, parity, nu- 



cleon numbers and so on. Then we restore all the symme- 
tries by applying symmetry projection (or in other words, 
quantum-number projection) operators onto symmetry 
broken mean-field wavefunction. Resultant quantum- 
number projected wavefunction is known to be able to 
give a better description. 

Here we consider strongly correlated electrons on a lat- 
tice, which have symmetries as total spin, total momen- 
tum, and so on. In general, explicit construction of sym- 
metry imposed wavefucntion is quite complicated. For 
instance, wavefunction with a definite total spin needs 
complicated spin coupling among a large number of elec- 
trons. However, the projection technique enables us 
to easily handle symmetry imposed wavefunction. This 
projection method is well harmonized with the recently 
proposed path-integral renormahzation group method 
(PIRG) 0,11 which has been quite a powerful tool and is 
free of the notorious minus-sign problem in investigating 
strongly correlated electron systems. In this method, the 
ground state is described explicitly by superposition of 
basis states, which often break symmetries possessed by 
the Hamiltonian when the numerically manageable num- 
ber of the basis states, L, is limited. By applying the 
projection operator to these basis states, we can exactly 
treat the symmetry and extract the state with a specified 
quantum number. We show that such a quantum-number 
projection technique can extensively widen applicabil- 
ity of the PIRG in the following points: (1) Precision 
of the numerical calculation is substantially improved. 
(2) The quantum number of the ground state is exactly 
determined. (3) The extended PIRG by the quantum- 
number projection can handle excited states and spectro- 
scopic properties in addition to the ground state. Such 
low-energy excitations correspond, in nuclear structure 
physics, to the yrast state 0, , which means the lowest 
energy state with specified quantum numbers (for nuclear 
structure, angular momentum). 

In Sec. II, we formulate the method of quantum- 
number projection with examples of spin, spin parity, 
electron momentum and lattice symmetry. In Sec. HI 
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and Sec. IV, we discuss an implementation of the 
quantum- number projection to the algorithm of the path- 
integral renormalization group (PIRG) method. Then 
several different ways of implementation are proposed 
in the order of increasing elaboration and accuracy. In 
Sec. Ill, we present algorithms and applications from 
this quantum-number projection technique applied after- 
wards to the obtained PIRG wavefunctions. Next, in Sec. 
IV, we show algorithms of the quantum-number projec- 
tion performed simultaneously with the PIRG procedure, 
by which the lowest energy state with the specified quan- 
tum number is more efficiently extracted. We show that 
the present methods applied in Sec. Ill and Sec. IV very 
efficiently improve the accuracy of the energy estimate. 
We show examples in the case of the Hubbard model. 
We also show how the excitation spectra are obtained. 
In the example of the Hubbard model with geometrical 
frustration effect, the present method enables to obtain 
the ground state as well as excitation spectra, which can- 
not be obtained in the existing methods. In Sec. V we 
summarize the results. 



II. QUANTUM-NUMBER PROJECTION 

In general, a basis state described by single Slater 
determinant does not often satisfy definite symmetry 
properties. Therefore, it can contain many components 
with unfixed quantum numbers, most of which are un- 
necessary for considering the specific eigenstate of con- 
sidered system. Here we consider a method to project out 
a component with a given quantum number from such a 
symmetry broken basis state. 

Projection operator C is usually defined as £^ — C If 
we act C onto wavefunction \il>), C lip) contains a compo- 
nent with the considered symmetry. By such quantum- 
number projected bases, the corresponding projected ma- 
trix elements are evaluated by (-01 C \^p), {^p\ HC \^) and 
(^1 OC for norm, Hamiltonian and other physical ob- 
servable matrix elements, respectively, where H is Hamil- 
tonian and O means a physical observable. Note that 
commutable property between observables and projec- 
tion operator and projection property — L simplify 
projected matrix elements. For the physical variables, 
we assume that O and L are commutable each other. In 
this section, we discuss the spin, momentum and lattice 
symmetries. 



A. Spin projection 

Quantum mechanically, finite object with a fixed shape 
must be rotated to recover the original symmetry. For 
nucleus, mean- field methods such as Hartree-Fock and 
Hartree-Fock-Bogoliubov approximations, give its opti- 
mum wavefunction. Though the rotational symmetry in 
the obtained wavefunction is broken, it directly relates 



the geometrical shape of nucleus. Restoration of rota- 
tional symmetry can be carried out by superposing ro- 
tated wavefunction. This superposition can be exactly 
carried out by angular momentum projection. Rotation 
in three dimensional space is specified by the Eulcr's an- 
gles and the restoration of the symmetry is usually de- 
scribed by the integration over the Euler's angles and 
weight of such superpositions is given by Wigner's D 
function. Angular momentum projection can be achieved 
by three-fold integration over Euler's angles as we will 
show later. Though this derivation is shown in nuclear 
structure textbook Q, in Appendix A, we discuss some 
properties of the projection operator. 

Here we first consider the spin degrees of freedom of 
electrons. Though the spin has no relation to any defi- 
nite shape, algebraic structure is the same. As the deriva- 
tion of angular momentum projection relies on the SU(2) 
structure, the same technique can be applied to electron's 
spin coupling. We consider to pick out the total-spin S 
component from a basis state described by a Slater deter- 
minant. As the Slater determinant has a definite number 
of up and down electrons {N^ and N\), z-projection of 

the spin is A'o — ^~ . This fact simplifies a projec- 
tion operator to a rather simple one. In nuclear structure 
physics, it corresponds to the case of angular momentum 
projection for axially symmetric shape. 

The spin projection operator has a form as 

2 9+1 t 

^ J dnDf,*^{m{^), (1) 

where — (a,/3,7) is Euler angle and Dffj^{Q) is 
Wigner's D function. Here M and K specify the z com- 
ponent of the total spin, Sz- As explained in the Ap- 
pendix A, Eg . (|A5(I . this projection operator operating as 
^MK IV') to a state \tp) filters out K component of |V') 
and generate a state which has Sz — M hy rotation. The 
rotation operator R{i}) is defined as 

= e*"^^e*'^'^''e*''^% (2) 

where 5*1, and Sz are y and z components of spin opera- 
tor, respectively. Wigner's D function is defined by this 
rotation operator as 

Df,i,{n) - {SM\R{n)\SK) = e^"^'^e^^^d£,^(/3), (3) 

where dfjj^i/i) = {SM \e'''^'^y\SK). By this projector, 
the spin projected state is written as 

LiiKm=LMNo\i^), (4) 

where Nq ^ {N^ - Ni)/2. Note that 1-0) has a definite Sz 
value. No, but e^^^" generates different Sz components. 
Therefore successive e^"^' selects finally needed Sz com- 
ponents. Although the Sz value is not unique and can 
have values in the range j^zj < 5 in the case of 5 / 0, 
this degree of freedom J| is eliminated by the following 
property of the spin projector; 

LmkLIpk' ^ Lfij(,Sss'SKM'- (5) 
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This relation can be easily proven by Ea. (jA5|) in Ap- 
pendix A. This relation shows that spin projection op- 
erator satisfies an extended projection property. As the 
PIRG basis states have a definite z-component of spin, 
the following relation is satisfied; 



tS tS 
^NaM^M'No 



(6) 



as the special case of Eq. (0 • Here we note that Lf^^ 
has a simpler form, which involves only one-dimensional 
integral, as 



^NoNo 



25 + 1 



d/3sin/34„^„(/3)e^''^='. (7) 



^NoNo suppressing S'^ 



In eq. jnj, we can take A^o as the M value. In this case, 
as the spin projection operator, we can use L^^^^ which 

satisfies usual projection property (Lff^^^j'^ = Lfj^-,pf„- 
Therefore in a later discussion, the spin projection oper- 
ator is simply denoted as — '^^ 
value. 

Because and H commute each other, 

{i^'\£^HC^\i;) = {il;'\H{C^f\iP) = {il;'\HC^\il;) is 
satisfied. Consequently, norm, Hamiltonian and 
other physical observable matrix elements between 
spin-projected basis of 1-0'} and [ip) are shown as 

^] 2S + 1 r f M 

H \ = — I d/3sin/34^^J/3)(V''| <^ H } 



O 



O 



(8) 



where rotated basis in spin space is defined as 

\^{/3))^e^PSy\i;). (9) 

Here we assume that O is a scalar operator for spin ro- 
tation and 5^ and O commutes 0. Note that, lip) is a 
direct product as lip) = \ tpi) \ tpi) , while its rotated one 
needs a larger representation space as the up and down 
components are mixed. 

For the case that the electron numbers of up and down 
spins are the same, the d function simply reduces to Leg- 
endre function Ps {cos P), 



4o(/3) =^s(cos/3). 



(10) 



Involved integral in eq. Q can now be efficiently evalu- 
ated by the Legendre-Gauss quadrature in practical nu- 
merical calculations. This quadrature needs less mesh 
points than those of trapezoidal formula. Typically, for 
S' = of the half- filled electron system in 6x6 and 12x12 
lattices, we needs 12 and 24 mesh points, respectively, for 
numerical convergence. As spin goes up, larger number 
of meshes is needed. 



B. Spin-parity projection 

We consider partial spin projection for the restricted 
case that the electron numbers of up and down spins are 



the same. Although it is not general, its scope is still 
wide. 

Now we consider the interchange between up and down 
spin components and define a parity for this interchange. 
We show that the parity classifies the even and odd total 
spins. Hereafter we call it spin parity. 

The parity operator may be defined as P = 
exp(— iTrS'y) = —iSy, where we obtain 

{SO\exp{~inSy)\SO) = 4„{7t) = PsIcostt) = 

(11) 

This reads that -I- parity wavefunction corresponds to 
even values for S and — parity wavefunction does to odd 
values. Therefore, this spin parity projection 



/:^± = (i±p)/2 



(12) 




yields to the classification between even and odd total- 
spin states. 

The spin-parity projected matrix elements are shown 

by 



(13) 



where with a — +1 and —1 takes [0) and P|0), 
respectively. 

If we take the spin projection operator, the spin-parity 
projection becomes redundant. However, in the case of 
multiple quantum-number projection operators, numeri- 
cal calculations inevitably become heavy. Since the whole 
spin projection is much more computer-time consuming, 
the spin-parity projection is an alternative way partic- 
ularly for the method of simultaneous quantum-number 
projection in each step of PIRG as proposed in Sec. IV. 



C. Momentum projection 

In systems with translational invariance, the conser- 
vation of momentum holds. However, a basis state is 
not necessarily an eigenstate of the momentum operator. 
By the projection technique, we restore the translational 
symmetry. We define the momentum projection operator 
as 



M ^ 

3 



(14) 



where A/" is the normalization, K is the momentum op- 
erator and Rj is a shift in a lattice specified by j. By 
applying this projection operator, we can calculate pro- 
jected matrix elements as 





\<t>{j)), (15) 
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where \(j> {j)) is a shifted wavefunction by the shift j. In 
an Lj. X Ly lattice, the momentum projection requires 
Lx X Ly larger computation efforts than those of unpro- 
jected one. 



D. Lattice symmetry projection 

In the Hubbard model on a two-dimensional lattice, 
there are several geometrical symmetries on a lattice as 
x-reflection, y-reflection and x-y interchange symmetries. 
Their symmetries can be classified by parity. By the 
associated parity operator P, we can define the corre- 
sponding projection operator as £ — similarly to 
the spin-parity projection. 



square lattice. We impose the periodic boundary condi- 
tion. 

We decompose exp[— ri/] into exp[— ri7] ^ 
[exp[— Ari/jj-] exp[— AriJj/.]]^ for small At, where 
T — MAt. When we use the Slater determinant as 
the basis functions, the operation of exp[— Ar_ffi^] 
to a Slater determinant simply transforms to an- 
other single Slater determinant. On the other hand, 
the operation of exp[— Arifj/.] can be performed by 
the Stratonovich-Hubbard transformation, where a 
single Slater determinant is transformed to a linear 
combination of two Slater determinants. 

One of numerical realizations of Eg. ljlGI) is PIRG 
method 0, 0. After the operation of exp[— ri7], the 
projected wavefunction can be given by an optimal form 
composed of L Slater determinants as 



III. QUANTUM-NUMBER PROJECTION TO 
THE PIRG STATES (PIRG+QP) 

A. Algorithm 

We briefly introduce the path integral renormalization 
group (PIRG) method, which has recently been proposed 
for solving strongly interacting electron systems 0, 13 . 
In general, the ground state \ipg) can be obtained by ap- 
plying the projector e~'^^ to an arbitrary state |(/>initiai) 
which is not orthogonal to the true ground state as 

\^Pg) = lim e-^^|Vi„itiai). (16) 

r — *oo 

In this paper, we consider the standard Hubbard model 
on a two-dimensional square lattice defined as 



H = Hi. 



where 



Hk — Ht + Hf 



H.c. 



and 



{kl)a 



Hut = U[n^^ - ]A ■ (mi - i 



(17) 



(18) 



(19) 



(20) 



(21) 



Here i, j represent lattice points and (cja) is a cre- 
ation (annihilation) operator of an electron with spin 
a on the i-th site. The summation over (ij) is for the 
nearest neighbor pairs and that over (kl) is for the next- 
nearest neighbor pairs on the 2D Hubbard model on the 



(22) 



a=l 



where Cq,'s are amplitudes of Operation of the 

ground-state projection can give optimal Cq,'s and |0q^^)'s 
for a given L. Its detailed algorithm and procedure are 
found in Ref. . 

By a finite number L, in most cases, it gives an over- 
estimate of the exact energy eigenvalue, since this wave- 
function satisfies the variational principle. Therefore, a 
relation between energy difference SE and energy vari- 
ance AE may be useful to extrapolate the energy into 
the true one. Here the energy difference is defined as 



SE = {H) - {H}g 
and the energy variance is defined as 



A^; 



H 



(23) 



(24) 



Here, {H)g represents the true ground-state energy. For 
we evaluate the energy and energy variance 
A£'(^), respectively. 

If is a good approximation of the true state, 

the energy difference SE^-^^ is proportional to the en- 
ergy variance AE^^\ Therefore extrapolating E^^^^ into 
^ by increasing L systematically, we can esti- 
mate accurate ground-state energy. 

Now we consider an implementation of the quantum- 
number projection to the state obtained by PIRG. The 
PIRG gives approximated wavefunction for a given L 

which is composed of L linear combinations of (f)'^^'^. 

One possibility to implement the quantum-number pro- 
jection is to project out as 



(25) 
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where /! is a quantum- number projection operator. We 

use the same amphtudes Cq, 's and the same bases (j)^a^ 

which the PIRG determines. On the other hand, this 
ampHtude c^'s can be easily reevaluated by diagonaliza- 
tion by using quantum-number projected bases, that is, 
we determine Cq's by solving the generalized eigenvalue 
problem as 



Ha 13^ — ^a(3^^ 



(26) 



where N^^ = {^fsl C\(l)a), H^^ ^ {M HC\(t>a). The lat- 
ter procedure gives a lower energy eigenvalue. By adding 
this procedure for the PIRG basis, we evaluate the pro- 
jected energies and energy variances, £^proj £^nd AE^j.^^ 
for each L. We can estimate accurate energy by extrap- 
olating the projected energy into zero variance. 

As a result of the application of this procedure, there 
appear two new aspects. One is that the energy estimate 
becomes more accurate. In general, correlation energy 
comes from dynamical and symmetrical origins. Original 
PIRG seeks for better basis states which gain both corre- 
lation energies in a compromised way. On the other hand, 
by the quantum-number projection operator, correlation 
energy originated in the symmetry is exactly evaluated. 
Consequently, the projected energy becomes much lower 
than the unprojected energy at a given L. If we use suffi- 
ciently large L, both values are the same and become the 
exact ground state energy. In practical problems, how- 
ever, we have to use finite number L and exact energy 
is estimated by extrapolation. Therefore, at the same L, 
better energy is useful for better estimation of the exact 
energy. 

The second point is that it enables the evaluation of 
excitation spectra. If we use projection technique, evalu- 
ation of excited states with different symmetry quantum 
numbers becomes easier. The PIRG basis states for L 
still have components of excitations which most likely 
belong to low-lying excited states. By projecting out the 
component with different quantum numbers from that of 
the desired one, we can evaluate such excited states. We 
note the lowest energy state with the specified quantum 
number (namely, the yrast state) is obtained. 



B. Numerical Results — Comparison to the exact 
results — 

We demonstrate how the method of quantum- number 
projection procedure applied to the PIRG wavefunction 
works by comparing with the exact results. 

First we consider the half-filled case on 4x4 lattice with 
U/t = 4.0. Its exact ground-state energy is -29.62185. 
The extrapolated energy of the PIRG is -29.488, when 
we use the data up to L = 320. We note that the 
auxiliary-field quantum Monte Carlo (QMC) method Q 
with rather large t ~ 20 — 30 also gives a similar value 
to that of the PIRG. There is some discrepancy between 
this energy and the exact one. This discrepancy comes 



from the remaining contribution from the higher-spin 
states contained in the projected wavefunction both in 
the PIRG and the QMC calculations. To obtain the real 
ground-state estimate, we need much larger t. Spin pro- 
jection can remove it very efficiently. In Fig. ^ we show 
spin projected energies of L=8, 16, 32, 64, 128, 256 and 
320 are plotted as a function of energy variances. The 
energy variance becomes smaller for larger L. In fact if 
the correct ground state is given, the variance becomes 
zero. As these energies are well scaled linearly as func- 
tions of the energy variance when the variance is small, 
the extrapolation to the zero variance works well. The ex- 
trapolated ground state energy is -29.611, which is quite 
close to the exact one. This result can also be compared 
with the variational Monte Carlo calculation with the 
Gutzwiller projection |^, which gives -29.47 jlflj. The 
SU(2) symmetric Monte Carlo calculation gives much 
better estimate with a reasonable value of r ~ 20. 
This is similar to the PIRG with the spin projection. In 
this sense, exact treatment of spin quantum number is 
crucial in obtaining the exact ground state in an efficient 
way in the present case. 

In Fig. n projected energies with S — 1 ^ 3 are also 
shown as functions of the energy variance. In addition to 
the ground state, excited states with 5=1^3 have a 
good linear scaling. Thus we can evaluate energies of the 
excited states with different spins by the present spin pro- 
jection technique. This fact shows an essential advantage 
of the PIRG combined with the quantum-number pro- 
jection technique, if one compares with the other type of 
numerical methods including the Monte Carlo methods. 
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FIG. 1: Extrapolation of the total energy to the zero energy 
variance for the spin projection for S = 0, 1, 2 and 3 in the 2D 
Hubbard model with 4 by 4 lattice and the periodic boundary 
condition. L is taken up to L = 320. The parameters are at 
t = l,t' = and U — 4. Exact energies with corresponding 
spin are shown by open diamonds. 

We investigate these extrapolations more closely. For 
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S = 1 and S = 2, extrapolated energies are very close to 
the exact ones, while for 5 = 3, the extrapolated energy 
is, to some extent, deviated from the exact one. As the 
PIRG is the projection to the ground state, the obtained 
wavefunction represents the ground state approximately. 
Therefore, as the total spin increases, amplitudes oi S ^ 
components in the PIRG wavefunction are expected to 
become smaller, because such high-energy component is 
already efficiently eliminated out by the PIRG projection 
process. Therefore extrapolated energy for higher-spin 
(for example 5 = 3 state) is worse than those of lower- 
spin (for example, S — and S — 1 states), because 
the higher energy states are almost missing in the PIRG 
states. Moreover, at a fixed L, the variance becomes 
larger as the spin goes up, which makes the extrapolation 
worse. This also indicates that the quality of projected 
wavefunctions becomes worse. We propose an improved 
algorithm to solve this difficulty for excited states in Sec. 
IV. 

Next we consider the spin-momentum projection. For 
the even or odd S, k — (0,0) or A; = (7r,7r) is consid- 
ered, respectively. In Fig.|21 we plot the spin-momentum 
projected energies as functions of energy variances. A 
remarkable difference between spin projection and the 
spin- momentum projection lies in the precision of energy. 
The extrapolated energy of the ground state is -29.62166 
at 5 = and k = (0,0). The accuracy is one order of 
magnitude better than the case of the spin projection 
only. As we show in Fig. 13 the spin-momentum pro- 
jected energy at L — 320 is -29.61650 while the energy 
with the spin projection only is -29.60228 for the same L. 
With the spin-momentum projection, at the same L, the 
energy becomes lower and extrapolated energy becomes 
closer to the exact one than that of the spin projection 
only. 

The higher spin state at S = 3 with the spin- 
momentum projection, to some extent, has a better ex- 
trapolated energies than the spin projection only, while 
there still remains a tendency that the extrapolation be- 
comes worse as the total spin goes up or the excitation 
energy increases. To overcome this defect, we have to 
consider the PIRG with projected bases, namely QP- 
PIRG method. We wiU show the efficiency of QP-PIRG 
in Sec. IV. 

We next study the half-filled system at 6 x 6 lattice 
with U/t = 4. In Fig. ^ we show the extrapolations 
of spin projected and spin-momentum and lattice pro- 
jected energies as functions of the energy variance. We 
take the PIRG wavefunctions for various choices of L 
up to 256. For the spin projection, we can get the low- 
est energy states (yrast states) oi S — 0, 1,2,3 from the 
PIRG wavefunction. On the other hand, for the spin- 
momentum-lattice projection, we further resolve them by 
their quantum number associated with the corresponding 
symmetries as 5 = 0, 2 with k = (0,0) and 5=1,3 with 
k = (7r,7r). Consequently variances of each L wavefunc- 
tion become smaller. Moreover the slopes of the linear ex- 
trapolation in the plot of the energy vs. variance asymp- 




0.01 

(<eV-<E>^) / <E>^ 

FIG. 2: Extrapolation of the total energy to the zero energy 
variance for the spin and momentum projections for { S — 0,2 
and k = (0,0) ) and ( 5 = 1, 3 and = (tt, tt) ) in the 2D 
Hubbard model with 4 by 4 lattice and the periodic boundary 
condition. L is taken up to L = 320. The parameters are at 
t = l,t' = and (7 = 4. Exact energies with the correspond- 
ing spin and momentum are shown by open diamonds. 



totically obtained at large L for the spin-momentum- 
lattice projection are smaller than the data with the spin 
projection only. The spin-only projection shows a slight 
underestimate of the ground-state energy after the ex- 
trapolation, which is ascribed to an insufficient number 
of L in this case. The extrapolated ground-state energy 
of spin-momentum-lattice projection is -66.8822. For the 
sake of comparison, the SU(2) symmetric auxiliary field 
Monte Carlo calculation under the constraint of the spin 
singlet gives -66.87 ± 0.05 Within the statistical 

error of the quantum Monte Carlo results, these two re- 
sults agree well each other as we see in Fig. 0] From the 
extrapolation, the ground state energy is inferred to have 
better accuracy than the Monte Carlo data. 



Next we consider the excitation energies. The spin 
projected and the spin-momentum projected approaches 
give similar values for low- lying states, although the pre- 
cision is better for the latter algorithm. Spin projected 
excitation energy of 5 = 1 and 5 = 2 state is 0.082 and 
0.249, respectively, while the spin-momentum-lattice pro- 
jection gives 0.081 and 0.238 for 5 = 1, A; = (tt, tt) and 
5 = 2, fc = (0,0), respectively. The accuracy appears to 
be similar for larger system sizes. 



Hubbard U=4, 4x4, energy extrapolation 
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FIG. 3: Detailed comparison of extrapolation of the total 
energy to the zero energy variance for the spin projection 
and spin-momentum projection for S = ground state in 
the 2D Ifubbard model with 4 by 4 lattice and the periodic 
boundary condition. The parameters are ai t = l,t' — and 
U = 4. Exact energy with corresponding spin and momentum 
is shown by open diamond. 



IV. PIRG WITH QUANTUM-NUMBER 
PROJECTED BASIS (QP-PIRG) 

A. Algorithm 

In the previous section, we considered the quantum- 
number projection after the PIRG wavefunction is ob- 
tained for the optimization of the ground state. To study 
the properties of excited states, we can further implement 
an improved algorithm of the quantum-number projec- 
tion in the PIRG method. That is to perform the PIRG 
procedure itself by using the quantum-number projected 
basis. 

In general, the ground-state projector e^'^^ to \ip) can 
be applied to lower the energy even within symmetry- 
imposed restricted space. When the Hamiltonian pre- 
serves some symmetry given by the projection £, that is, 
when C and H are commutable, the lowest-energy state 
of the specified quantum number, 1-0), can, in principle, 
be calculated from 

IV') = lim e-"^/: I^Ainitiai) . (27) 

r — >oc 

By introducing the Stratonovich-Hubbard transforma- 
tion, however, a partial sum over the Stratonovich aux- 
iliary variable destroys the symmetry. Therefore, if one 
wishes all the time to keep the symmetry of the state 
with the specified quantum number, in an elementary 
PIRG procedure of the projection exp(— AriJ)|-0), we 
need to perform the quantum-number projection every- 
time as £exp(— Ar_ff)|?/') to restore the symmetry. This 
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FIG. 4: (color) Detailed comparison of extrapolation of the to- 
tal energy to the zero energy variance for the spin projection 
(blue open circles), spin-momentum projection (filled black 
circles) and spin-momentum-lattice projection (orange open 
circles) for S — 0, 1,2 and 3 in the 2D Hubbard model with 
6 by 6 lattice and the periodic boundary condition. The pa- 
rameters are sd, t — l,t' = and [7 = 4. Quantum Monte 
Carlo energy for the ground state is shown by open diamond 
with error bar. 



is a much more efficient way of obtaining the lowest- 
energy state with the specified quantum number than 
the PIRG-f QP method discussed in Sec. III. 

We here explain the algorithm more precisely in 
the case of the Hubbard model defined by Eq. H17|l . 
The basic procedure is then summarized as re- 
peated operations of C together with the operation of 
exp(— Arif). Namely, lim e~'^^ C\^pinitiai) is replaced 

r — >oc 

with hm [/:e-^^«^ ^^ Ce-^^""^]^ IV-initiai) by keeping 

r — >oo 

At small. Here the operation of Q-^'^^ui contains the 
Stratonovich-Hubbard transformation. A partial and op- 
timized sum of the Stratonovich-Hubbard auxiliary vari- 
able constitutes the truncation of basis to keep the num- 
ber of basis, while it destroys the symmetry. This al- 
gorithm allows the restoration of the required symme- 
try by the operations of C at each step of the trunca- 
tion. This is the best way of the optimization to obtain 
the lowest energy state which has the required symmetry 
(namely, the yrast states). In each step of the operation 
of exp[— Ar_ffif] or exp[— Ar_ffc/J, we employ the trun- 
cated basis which gives the lower energy for the states 
£exp[-Ar_ffK]|V') or £ exp[— ArTJ^/J IV"). We cah this 
algorithm of simultaneous PIRG and quantum-number 
projection, Quantum-number Projected PIRG (QP- 
PIRG). To differentiate from QP-PIRG, the quantum- 
number projection procedure using the original PIRG 
result explained in Sec. HI is called PIRG+QP. 

In principle, any quantum-number projection operator 
can be used in the PIRG. However, in practical appli- 



cations described later, we take a set of multiple pro- 
jections, namely spin-parity projection and momentum 

projection operators, C^^C'^. Ideally, all the quantum- 
number projection operators should be applied, while it 
rapidly increases numerical computation time. In the 
present paper, as we study the full momentum disper- 
sion, we employ the momentum projection operator. Al- 
though the spin projection is important, the spin rotation 
in spin space mixes the up and down spin components 
and we need the twice as large space as the original one 
for the Green function in the PIRG procedure, which 
makes the PIRG computation heavy. Then for the mul- 
tiple projection of QP-PIRG, to save the computation 
time, we propose, for a practical use, a combination of 
the momentum and the spin-parity projection instead of 
the full spin-momentum projection. By this approach, 
the PIRG wavefunction does not have a good spin quan- 
tum number. To restore the spin symmetry perfectly, 
after the QP-PIRG procedure above, we again perform 
the full spin projection afterwards. Namely, to obtain a 

final result, C'^ Ciatuce is applied. This constitutes the 
full procedure of QP-PIRG. 

At each quantum-number projection, the integrations 
or summation such as those in Eqs.® and (|15|l . can be 
very efficiently parallelized in actual computations if par- 
allel processors are available. In each process, we store 
the Green function Q^^"'^^/^^ = (^/;^|c|cj£|V'/3), while the 
update of the Green function after the operation of each 
g-Arffc/i jg -ftfritten as 

4^-^'^^)^(^„|clc,£|^^), (28) 

or 

G\f'^^^^^{^yjclc,C\i.,), (29) 

where 

= ^ iexp[2aa(n,T - n^) - AtU/2Up), (30) 

with a being the Stratonovich auxiliary variable and 
a = tanh^^ tanh(-^j^). When one term of the sum 
over a is taken in the truncation process, the updated 
Green function is efficiently calculated from the old Green 
function df'^'^'^"'^ in the same way as Eq.(3.10)-(3.14) 
in Ref. 8]. 

B. Numerical Results of QP-PIRG 

1. Results for 6 by 6 lattice 

Now we show numerical results of QP-PIRG. We first 
show the case of 6 by 6 lattice at [/ = 4 and t = l,t' ^ 0. 
In Fig. |S1 we show the extrapolation of QP-PIRG re- 
sult by green open circles by using the projection up to 
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FIG. 5: (color) Extrapolations of the energy to the zero en- 
ergy variance by using PIRG-I-QP (filled blue circles and dot- 
ted lines) and QP-PIRG (open green circles and solid lines) 
for the 2D Hubbard model with 6 by 6 lattice and the periodic 
boundary condition. The parameters are at f = = and 
U = 4. The ground-state energy of Monte Carlo calculation 
is also shown by open brown diamond at zero variance with 
the statistical error bar (-66.8664 ± 0.0504). The red symbols 
with red solid lines are derived from the largest L wavefunc- 
tion (L = 140) of QP-PIRG, where we choose partial La basis 
functions which have the largest weights among L bases. The 
plots are obtained with increasing La up to L = 140. 



L = 140. As we discussed, the QP-PIRG with quantum- 
number projected bases seeks for optimum yrast states 
concerning the considered symmetry in every PIRG pro- 
cess. In this calculation, we took spin-parity and mo- 
mentum projection operators. For S* = and k — (0, 0) 
state, we use £•^+£'==(0.0). As the obtained wavefunc- 
tion still contains S = 2,4,.. components, we apply 
£S=0£fc=(o,o)£^^^^.^^ projection operators afterwards for 
final results. 

This QP-PIRG can generate a better wavefunction 
than the PIRG+QP state as we see in the compari- 
son with blue closed circles. Here we show results of 
the PIRG-I-QP state obtained after spin-momentum pro- 
jection. In fact, for 5 = with k = (0, 0) state, in 
the PIRG-I-QP result even at L = 256, the energy is - 
66.5765, while the same energy can be given at L ~ 15 
by the QP-PIRG. This means that for the ground state, 
basis states are more elaborately selected by the QP- 
PIRG. Thus, the quantum-number projection simultane- 
ously with the PIRG provides an efficient way of obtain- 
ing better wavefunctions. The extrapolated ground-state 
energy is -66.879 which is well within the statistical er- 
ror of the previously cited Monte Carlo energy. In fact, 
from the extrapolation procedure in Fig. \^ the accuracy 
of the QP-PIRG seems to have more than 4 digits and is 
higher than the accuracy of the presently referred quan- 
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turn Monte Carlo result , since the energy at L = 140 
is already lower than the upper bound of the Monte Carlo 
estimate. Namely, the QP-PIRG result seems to give the 
highest accuracy among these comparisons. 

In addition, we have also shown in Fig. an alter- 
native way of the extrapolation. The red symbols are 
derived from the largest L state after QP-PIRG, where 
L = 140 in this case. This state is represented by L basis 
functions. After ordering these basis functions from the 
largest weight in the linear combination, we may trun- 
cate the basis functions by taking only the La states from 
that with the largest weight. By using these truncated 
functions with different Laj we have plotted the energy 
and variance of these truncated states. This gives very 
close estimate to the QP-PIRG result shown above as 
the open green circles. A small difference between this 
procedure and the original QP-PIRG is seen at larger 
variance. This may be due to the fact that at small L, 
the present truncation at small La does not necessarily 
give the lowest energy state with La- Another possible 
origin is that the iteration of the present QP-PIRG is not 
sufficient in reaching the lowest energy state under the 
constraint of each L. In any case, the linearity of the plot 
in the plane of the energy and the variance is well satis- 
fied in both cases, particularly for the latter procedure, 
and the asymptotic slopes at large L look the same. 



problem at this parameter value. We compare our results 
with the exact one. 
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FIG. 6: (color) The energy dispersion of the 2D half-filled 
Hubbard model at [/ = 5.7, t = l,t' = 0.5 for S = states. 
The system size is 4 by 4 with the periodic boundary condi- 
tion. The comparison with the exact results (black crosses) 
shows that the QP-PIRG (red circles) works excellently well 
for the ground state as well as the dispersion even when the 
geometrical frustration efi^ect is large. 



2. Results with next-nearest neighbor transfer 



In the previous section, we consider the standard Hub- 
bard model with t' = 0. Conventional quantum Monte 
Carlo calculation could be performed to investigate such 
ground state properties. The PIRG is an alternative 
method in this respect while it and its extension have 
an advantage in investigating the excitation spectra. Es- 
pecially, quantum-number projection enables us to han- 
dle yrast states with the same effort as the ground state. 
However, it is expected that the 2D Hubbard model with 
t' — Q has an antiferromagnetic long-ranged order in 
the thermodynamic limit and has a simple low-energy 
structure. To test the efficiency of our algorithm in a 
more severe condition, we investigate the extended Hub- 
bard model by including the next-nearest neighbor trans- 
fer, which causes the geometrical frustration effect. The 
quantum Monte Carlo method is known to have a severe 
difficulty when t' becomes large. 

Recently by using the PIRG method, the non-magnetic 
insulator (NMI) phase has been found near the Mott 
transition for relatively large t' fl5|. This phase can not 
be investigated by the Monte Carlo methods due to se- 
vere minus sign problems. Therefore, the PIRG is so far 
the only technique suited for this study. Here we ex- 
plore how the present quantum-number projection tech- 
nique improves the precision of the PIRG in such a study. 
Here we consider the half-filled system on 4 by 4 lattice 
with U /t = 5.7 and t' = 0.5. Monte Carlo method does 
not give us convergent results because of the minus sign 
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FIG. 7: (color) The energy dispersion of the 2D half-filled 
Hubbard model ai U = 5.7, t — l,t' — 0.5 for S — 1 states. 
The system size is 4 by 4 with the periodic boundary condi- 
tion. The comparison shows that the QP-PIRG (red crosses) 
works excellently well even for the spin excitations. 

In Figs.|n|and[3 we show comparison of the dispersions 
obtained by the QP-PIRG with the exact diagonalization 
result. This system has the ground state at 5* = and 
k = (tt, 0). The 5 = with k = (tt, tt) state severely com- 
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petes with this ground state. The lowest-energy S = 1 
state has k = (tt, 0). This energy is very close to those 
of doubly degenerate k = (0,0) states and the second 
lowest state with 5 = 1 and k = (tt, 0). The comparison 
of QP-PIRG (red circles) and the exact diagonalization 
results (black crosses) in Figs. El and [3 indicates excel- 
lent agreement. In general, the errors are less than 0.01, 
which means the accuracy higher than the 4 digits. 

V. SUMMARY 

We have presented the quantum-number projection 
technique and its implementation to the PIRG method, 
both of which works well irrespective of the details of 
the considered system. The quantum-number projection 
method can pick up a component with required symme- 
tries from symmetry broken wavefunctions (i.e., mean 
field wavefunction and so on). In the Hubbard- type 
model, the symmetries have a significant role in the low- 
energy states. In particular, spin, momentum and lattice 
symmetries play specially important roles in determining 
the low-energy states. Restoration of the spin symmetry 
can be carried out by taking a spin projection operator, 
which is the same technique as the angular momentum 
projection in nuclear structure physics. Spin rotation is 
performed in the spin space and the spin projection is 
represented by one dimensional integral for the rotation. 
The momentum projection is simply given from the su- 
perposition of spatially translated basis functions. We 
have also considered geometrical symmetry on a lattice 
for projections such as the inversion and rotation sym- 
metries. 

Quantum-number projection operator C is represented 
by the sum of exponential of one-body operator. In the 
PIRG, the wavefunction is expressed by a linear combi- 
nation of L basis states, while the symmetries are not re- 
tained in each basis state in general. Then the quantum- 
number projection is efficiently introduced for each ba- 
sis state. In the present paper, we have introduced two 
ways of implementing quantum-number projection into 
the PIRG. One way is to carry out quantum- number pro- 
jection afterwards for the already obtained PIRG wave- 
function (PIRG-I-QP). The ground state is efficiently ex- 
tracted by specifying the quantum number with higher 
accuracy than the PIRG only. Although the PIRG does 
not efficiently pick up the excited states, we can obtain 
several low-lying excited states with various symmetries 
from the PIRG wavefucntion, if a small portion of the ex- 
cited states still remain after the PIRG procedure. Other 
is to carry out the PIRG by using quantum-number pro- 
jected basis states (QP-PIRG). By this extended PIRG, 
we can precisely evaluate excitation spectra. Although 
QP-PIRG requires more computation time, the accuracy 
of the ground state is more improved than PIRG+QP, 
particularly for the excitation spectra. 

In numerical calculations, quantum-number projection 
can be performed exactly. Moreover, as quantum- number 



projection operators C are commutable with Hamilto- 
nian if, the relation CHC = HC simplifies numerical 
calculations. As examples, the accuracy and efficiency 
of the algorithm has been tested for the standard Hub- 
bard model on two-dimensional square lattice as well as 
for the 2D Hubbard model with nonzero next nearest 
neighbor transfer, where geometrical frustration effects 
are large. We have shown that the quantum-number 
projection implemented to the PIRG excellently works. 
More concretely, the spin projection and spin-momentum 
projection by PIRG-f QP greatly improve the accuracy of 
energy. QP-PIRG further improves the accuracy of the 
extrapolated energy. This algorithm also enables accu- 
rate calculations of low-lying excitation spectra with dif- 
ferent quantum numbers from those of the ground state. 
The energy dispersions of the specified total spin have 
been shown to give highly accurate results, particularly 
by using the QP-PIRG method. This accuracy does not 
depend on the details of the lattice structure and the 
dimensionality. In our examples the accuracy becomes 
higher or comparable to 4 digits. 

When the system size increases in the 2D Hubbard 
model, we do not have a relevant clue to judge the accu- 
racy of the calculation by the present algorithm. On the 
half-filled case, however, we can compare the results with 
the quantum Monte Carlo results and the agreement is 
satisfactory. To reach the same accuracy, it seems to be 
necessary to increase the number of the basis functions 
L gradually with the increase of the system size. 
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APPENDIX A 

In this appendix, we discuss some properties of the 
spin projection operator. 

We expand by complete set \SMa) regarding to 
spin quantum number, as 

IV") = csMo.\SMa), (Al) 

SMa 

where csMa = {SMa\tp) and a denotes other quantum 
numbers. Operation of rotational operator R{U) to |^/;) 
results in 

Rm\^) = J2 ''SMo.RmiSMa) 

SMa 

= csMaDi^t{n)\SKa), (A2) 

SKMa 
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where we use eq.Q. By this relation, projection onto 
is represented by 



Therefore, Lfj^^ projects out \SM) component from 
By Ea. (|A3() . projection operator is represented as 



25+1 



87r2 

Y^\SMa){SKa\il^) 



(A3) 



where we use the foUowing relation as 



'^\SMa){SKa\ 



(A5) 



'25 + 1 
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